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ABSTRACT 

The first stars forming in minihaloes at redshifts greater than 20 may have been very 
massive and could have left behind massive black hole (MBH) remnants. In previous 
papers we investigated the hierarchical merging of these 'seed' MBHs and their as- 
sociated haloes, using a semi-analytical approach consisting of a hierarchical merger 
tree algorithm and explicit prescriptions for the dynamics of merged substructure in- 
side a larger host halo following a merger. We also estimated accretion luminosities 
for these MBHs and found them to be consistent with observations of ultra-luminous 
iy) , X-ray point sources. Here we compute the strength of gravitational wave events as 

Q\ ■ MBHs merge to form the more massive black holes that we predict reside in galaxy 

haloes today. If MBHs merge efficiently, we predict that as many as 10 4 -10 5 events 
C*~) . per year may fall within the sensitivity limits of the proposed LISA gravitational wave 

observatory. The collapse of the first massive stars to form MBHs may also be accom- 
^ (— i panied by gamma-ray bursts (GRBs). If this is the case and if GRBs are observable 

O out to the redshifts of first star formation, we predict that about 10 5 -10 6 GRBs per 

I 1 year could be detected. As merging MBH binaries reach their last stable orbits before 

final coalescence a fraction of the gravitational wave energy may be released as a pulse 
■ of gamma rays (for instance, through interaction with material enveloping a merging 

^ | MBH binary). This fraction has to be larger than about 10~ 2 for MBH mergers to 

account for some beamed GRBs, and greater than 10 -6 for the gamma rays to be 
detectable out to cosmological distances with upcoming GRB detector missions. 

' Key words: galaxies: formation - galaxies: haloes - galaxies: nuclei - cosmology: 

. theory 



1 INTRODUCTION 

There is increasing speculation that the first generation of 
stars in the universe may have been extremely massive, and 
that some of these objects could have collapsed directly to 
massive black holes (MBHs) at the end of a brief stellar 
lifetime. If these objects do form in the early universe, then 
they may provide a unique opportunity to study primordial 
star formation, through their gravitational wave or gamma- 
ray emission. 

Recent semi-analyti c (|Hutc hings e t al. 2002] 
IFuller fc Couchman 2007)1 ITegmark et al. 1997) and 

numerical studies (jBromm^T^ojjpTT^r^arson 2002| 

|Abel, Bryan fc Norman 2000} suggest that the first 
stars in the universe formed inside dense baryonic cores, 
as they cooled and collapsed within dark matter haloes at 
very high redshift. For a standard ACDM cosmology, these 
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minihaloes are estimated to have had masses in the range 
Mmin ~ 10 5 -10 6 Mq, and to have collapsed at redshifts 
^collapse ~ 20-30 or higher. Since these first star-forming 
clouds contained essentially no metals, gas cooling would 
proceed much more slowly in these systems than in present- 
day molecular clouds. As a result, they may have collapsed 
smoothly and without fragmentation, producing dense cores 
much more massive than the proto-stellar cores observed 
in star-formation regions today. Assuming nuclei within 
these cores accreted the surrounding material efficiently, 
the result would be a first generation of protostars with 
masses as great as 10 3 Mq ( |Bromm, Coppi fc Larson 2002| 
IQmukai fc Palla 20011 . 

As yet, nothing definite is known about the subsequent 
evolution of such objects. However, their large masses would 
probably result in many of them ending up as MBHs with 
little intervening mass loss - for systems in this mass range, 
gravity is so strong that there is no ejection of material from 
the system in a final supernova bounce ( |Heger et al. 200*21 1. 
This high-redshift population of MBHs is interesting as a 
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Table 1. The mass of seed MBHs, and the height of the peaks in 
the initial density field in which they formed, for the four models 
considered. The collapse redshift z co n apae is the epoch when peaks 
of height v pk first cross the cooling threshold (see paper I). 



Model 




peak height f p j. 


■^collapse 


A 


260 M 


3.0 


24.6 


13 


1300 M 


3.0 


24.6 


C 


260 M Q 


2.5 


19.8 


D 


260 Mq 


3.5 


29.4 



source of seeds for the formation of the supermassive black 
holes (SMBHs) seen at the centres of galaxies at the present- 
day. Furthermore, the initial collapse of these objects or 
mergers between them might be detectable through the 
resulting gamma-ray or gravitational wave emission. This 
would provide a direct test of the physics of primordial star 
formation and the high-redshift universe. 

In a previous paper ( Islam , Taylor fc Silk 2003, herafter" 
we used a semi-analytical model of galaxy formation to 
follow the evolution of MBHs as they merge together 
hierarchically, along with their associated haloes. The code 
we used combines a Monte-Carlo algorithm to generate 
halo merger trees with analytical descriptions for the main 
dynamical processes - dynamical friction, tidal stripping, 
and tidal heating - that determine the evolution of merged 
remnants within a galaxy halo. We introduce seeds into the 
code by assuming that in each minihalo forming before a 
redshift z C oiia P so, a single MBH forms as the end result of 
primordial star formation. 

For our computations, we considered four different sets 
of values for the parameters v p k and M, iSee d, which fix the 
abundance and the mass of the seeds respectively. These 
values are summarised in tabled Our choice of a seed MBH 
mass of 260 Mq is motivated by the result of Heger (2002) 
that massive stars above this mass will not experience a 
supernova at the end of their lives, but instead will collapse 
directly to a MBH of essentially the same mass. 

In paper I we investigated the abundance of MBHs 
in present-day galaxies as a result of this process of hi- 
erarchical merging and dynamical evolution. For these 
MBHs we estimated the X-ray accretion luminosities us- 
ing two different accretion models and found that they 
could well account for many of the ultra-luminous off- 
centre X-ray point sources that are observed in local 
galaxies( |Islam, Taylor fc Silk 2003, herafter paper II| | . 

In this paper we investigate how our model could be 
tested at high redshifts. In particular, we predict the number 
of gravitational wave events arising from MBH mergers that 
should be detectable with the upcoming LISA mission. Con- 
versely, observations of gravitational wave events could then 
be used to constrain the merging history of MBHs and pro- 
vide limits on the abundance of these objects at cosmolog- 
ical distances and to very high redshifts (see also Haehnelt 
(1994) , Menou, Haiman & Narayanan (2001) ). If the col- 
lapse of the first massive stars to form MBHs is accompanied 
by gamma-ray bursts (GRBs), then observations of GRBs 
may allow us to directly probe the epoch of first star forma- 
tion at redshifts larger than 20. Here we estimate the number 



of GRB events that could be observed per year if GRBs are 
detectable out to the redshifts of first star formation. 

The structure of the paper is as follows. In section [5] 
we compute the number of MBH mergers and the number 
of associated gravitational wave events. Section [3] consid- 
ers the expected number of GRBs from collapsing massive 
population III stars, as well as the possibility of gamma- 
ray emission in the wake of MBH mergers. We discuss and 
summarise our findings in section [I] 



2 GRAVITATIONAL WAVES FROM 
SMBH-MBH MERGERS 

There are essentially two types of gravitational waves that 
are emitted at different stages in the evolution of the sce- 
nario we are considering. The first of these are 'burst' signals 
associated with the collapse of massive population III stars 
into MBHs. In our model, these are events occurring at very 
high redshifts > 20 for a cosmological abundance of massive 
population III stars, so we expect a stochastic gravitational 
iwtve background to be generated llSchneider et al. 20001 
|de Araujo, Miranda fc Aguiar 2002^ . Although statistical 
limits can be imposed, the particular number and strength 
of burst signals created in this way does depend on the ini- 
tial number and masses of collapsing stars and the collapse 
mechanism. 

In contrast periodic signals, albeit with changing fre- 
quency, are emitted during the subsequent merger of MBHs 
to form more massive BHs. Compared to the burst signals 
from collapsing population III stars, the mergers of MBHs 
produce stronger waves due the higher total mass of the 
gravitationally interacting system. As a result gravity waves 
could be detected for individual merger events. The signals 
from merger events can be further subdivided into signals 
associated with the inspiral phase, actual coalescence and 
subsequent ringdown (Flanag an fc Hughes 1 998). Of these, 
the inspiral phase is the longest and thus most likely to de- 
tect. In the following we will be exclusively concerned with 
the gravitational waves emitted during this phase. These 
are also referred to as 'chirp' signals, due to the fact that 
during inspiral the distance between the BHs in the merg- 
ing binary decreases continuously as energy is lost through 
gravity waves, resulting in the frequency of emitted waves 
rapidly sweeping upwards. 



2.1 Gravitational Wave Amplitude 

For any measurement of gravitational waves the most im- 
portant pieces of information are the frequency and some 
measure of their amplitude, both of which will be functions 
of time. Because of the spin-2 nature of gravity, the first non- 
zero and largest term of the tensor describing the gravity 
wave induced perturbation of space-time is the quadrupole 
moment. This can be described by an overall dimension- 
less strain amplitude h s . For periodic signals such as those 
emitted by a binary black hole system, this is given by 
ornc 1989) 
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where f gw is the frequency of gravitational waves, which in 
this case is just twice the orbital frequency of the binary; 
M12 = Mi + M 2 and fii 2 = M1M2/M, G is Newton's con- 
stant and c the speed of light. D(z), the distance to the 
binary, can be obtained by integrating the differential dis- 
tance redshift relation 



dz' 



(2) 



The effective strain amplitude measured at some charac- 
teristic frequency f c is further enhanced proportionally to 
the square root of the number of periods, n, that the bi- 
nary emits at or near that frequency. For an inspiralling 
binary system whose frequency sweeps up, or 'chirps', in 



frequency (cf. below), n oc /j,^ M^'" f c °'" and the 
characteristic strain amplitude, h c is thus HThorne T9891 
INakamura et al. 19971 
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where f c is the frequency of the wave signal, r the distance 
from the merging binary and the 'chirp' mass M c hir P ac- 
counts for the masses Mi and M 2 of the binary constituents 
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Particularly at late times, most mergers in our model 
are between a massive central MBH and much lighter seed 
MBHs, so we can consider the problem in the limit that one 
of the merging binary constituents is much larger than the 
other. If we assign Mi and M 2 to the masses of the central 
and inspiralling MBH respectively, Mi S> M 2 and 
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where x = M 2 /Mi denotes the mass ratio between the two 
masses. 

The largest frequency at which a gravitational wave 
source can emit is determined by its light radius which in 
turn is dependent on its mass, and so in this case f ma x ~ 
10 Hz(lO 3 M0/Mi). A lower limit on the characteristic wave 
amplitude is then 
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2.2 Frequency of gravitational Waves 

The proposed LISA (Laser Interferometer Space Antenna) 
gravitational wave observatory will be the only instrument 
capable of detecting gravity waves in the typical frequency 
range generated in the inspiral phase of a merging binary 
of massive black holes. A rough estimate of the maximum 
frequency of waves emitted by the merging binary is 



fn 



10 Hz 



1O 3 M 
Mi 



(7) 



where Mi is the more massive of the binary constituents. 

For a binary on circular orbit, the frequency of the chirp 
signal can be described analytically. If the orbital radius is oo 



at time to = 0, then the time for the binary constituents to 
spiral into each other is (Mi sner, Thorne fc Wheeler 1973} 



to = 0.665 
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The frequency of the waves emitted at a later time t is 



f gw = 3.11 x 10 
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Including the redshift z, the waves are then detected by 
LISA with a frequency 
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For the case that Mi ^> M 2 , to and f gw are approxi- 
mately given by 
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2.3 MBH Merger Efficiency 
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We have outlined above the main characteristics of the gravi- 
tational wave signals received from individual merger events. 
The formulae given estimate the amplitude and frequency 
averaged over all inclinations of the merger orbital plane. 

For a merger event to actually occur the MBHs have 
to come within a distance where gravitational waves can 
efficiently reduce the orbital energy and allow the MBHs to 
spiral together and coalesce. However, this only happens at 
MBH separations that are very much smaller than the ones 
we are considering here (of order 10 - 100 pc). 

For a MBH binary on a circular orbit to coalesce 
through gravitational wave emission within a Hubble time, 
the MBH separation, a, needs to satisfy IPeters 1964t 



a < a a 



256 thubbie G z Ml ^ 1/4 
5 c 5 



.. o°(AA l/ V mi V- 
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At distances larger than this the common environment 
of the MBHs determines whether and how quickly the MBHs 
will be driven towards each other. If two MBHs come within 
a distance 

G M 2 



a < ah 



4 al 



(14) 
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pc 



they form a 'hard' binary. Interactions with stars in the bi- 
nary's vicinity carry away energy from the binary with the 
result that the binary separation shrinks - the binary is said 
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to 'harden'. The time scale on which this happens is given 
by 



G p* a H h 



(15) 



where Hh is the dimensionless hardening rate and is of 
the order Hh ~ 15 (- Quintan 1996| l. In a fixed isothermal 
stellar background with a number density of stars, p+ — 
cr 2 /(2 irGr 2 ) the hardening timescale starting at a 
becomes 



a h 



t h = 



7T G M 2 
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yr 



As it stands this implies that the hardening timescale is 
proportional to the lighter MBH in the binary system. This 
makes sense: as the mass of M 2 increases the fraction of 
energy lost to any interacting star decreases. A light MBH 
binary therefore needs less stellar interactions to harden. 
However, this is only strictly valid if the stellar background 
is fixed. That is, we have ignored any depletion of the stel- 
lar density in the binary environment due to the interactions 
and resulting ejections of stars from the binary (see for ex- 
ample Milosavljevic & Merritt (2001) ). 

The hard binary stage can represent a bottleneck on 
the way to MBH mergers - particularly for very mas- 
sive (S)MBHs, as dynamical friction is no longer signifi- 
cant but gravitational radiation not yet effective enough in 
further reducing the separation between the MBHs. How- 
ever, both, the dynamical interaction with stars in the MBH 
binary's vicinity as well as gas infall are likely to have 
the potential to quickly reduce the distance between the 
MBHs fBarnes fc Hernquist 1996| INaab fc Burkert 20011 
IBarnes 20021 |Armitage fc Natarajan 2002| lEscala et al. I 
We do not quantify this effect, but note that the role of gas 
infall in this context is probably more important at early 
times when MBHs encounter each other more often in the 
wake of major mergers of two host haloes 1 . 



2. 3. 1 Dynamical friction at the centre of haloes 

At still larger distances, dynamical friction is the main mech- 
anism by which MBHs are drawn to the centre. The dynam- 
ical friction force is given by the Chandrasekhar formula 
| |Binney fc Tremaine 1987| 



dt 



= -167r 2 lnA c G 2 m(M sat + rri 



f(v')v' 2 dv 



(17) 

where v sa t is the satellite velocity with respect to the host 
centre, M sa t the satellite mass, m is the mass of individual 
host halo particles, that have a distribution function f(v), 
and In A c = ln(& max /6 m i„) is the Coulomb logarithm, that 
is the ratio of maximum and minimum impact parameters, 



1 At very high redshifts z > 10 the mass ratios M2/M1 between 
MBHs are closer to 1 on average as central SMBHs have not yet 
grown to the massive sizes observed today. The more similar mass 
of MBHs implies that they encounter each other in the wake of 
the merger of two haloes that are also of similar mass. 




Figure 1. Ratio of dynamical friction time to orbital period 
plotted against ratio of host mass interior to satellite orbit to 
satellite mass. From top to bottom the curves are for values of 
the circularity parameter e = 1.0, 0.8, 0.6, 0.4, 0.2. 



satellite system. 

This effect has been modelled in the numerical scheme 
using the semi-analytical dynamical friction timescale given 
by Colpi, Mayer & Governato (1999; see also Taylor & Babul 
2004) 



r df = fc- 



(M hoa t/M aat ) 

ln(Mhost/M sat ) 



^ (18) 

I TV 



where Mhost and M sa t denote the mass of the host and satel- 
lite haloes respectively and P„; r is the circular orbital period 
at the virial radius of the host, and fc is a numerical con- 
stant. Colpi et al. (1999) find fc = 1.2e ~ 3.26 for massive 
satellites in an isothermal potential with a constant-density 
core. We will use the value k = 2.4, found to match the in- 
fall times for massive satellites in a cuspy potential in the 
model of Taylor and Babul (2004). The circularity parame- 
ter e = J(E) I Jdrc(E) is the ratio of the angular momentum 
of the actual satellite orbit and that of a circular orbit with 
the same energy. 

In the following we assume that this timescale domi- 
nates both the timescale for gravitational wave induced co- 
alescence and that for overcoming the hard binary stage, 
which in principle can be very short due to the effects of gas 
and stellar dynamical processes. 

In order to compute the time when a MBH merges with 
the central (S)MBH in the host we first establish the time 
when a satellite/MBH system sinks to an orbit with a ra- 
dius less than the radius of the infall region. This is 1 per 
cent of the host virial radius at the time. At this point the 
satellite/MBH system in question has so far been considered 
as having 'fallen to the centre'; following their dynamics to 
smaller radii was too expensive computationally. 

To track infall more accurately, we add to this time a 
dynamical friction infall time, using eg. 11811 and replacing 
the total host mass interior to the virial radius and the or- 
bital period at the virial radius with the respective quantities 
at the outer radius of the infall region. This is also shown 
in figure for a range of satellite to host mass ratios. The 



MBH remnants of the first stars III: observational signatures from redshift evolution 5 




10 3 r 



io 7 10 s 10' io 10 10" io 12 

M i,™ [ M «u,l 

Figure 2. The figure shows the total mass of the host halo versus 
the mass of the central (S)MBH. For the (S)MBH mass we have 
taken account of the slightly altered growth history due to the 
addition of a dynamical friction infall time after satellite MBH 
systems have crossed the infall radius. 



circular parameter e is computed from the last orbit of the 
satellite system before it crossed into the infall region. The 
result is a more accurate estimate of the time when the ac- 
tual merger occurs. As a result of the correction to the infall 
time we also inevitably change the order in which the MBHs 
merge with the centre. This in turn leads to a change in the 
growth history of the central (S)MBH and also to a slightly 
lower final mass of the SMBH, since there is now a small 
fraction of satellite MBH systems that do not merge with 
the centre within a Hubble time. Figure [5] shows the host 
halo mass versus the resulting mass of the central (S)MBH 
for final halo masses of 1.6 x IO 10,11,12 M and their lighter 
precursor haloes at higher redshifts 2 . The data points from 
precursors of different final halo masses are recognisable by 
their accumulation near those final halo masses. Neverthe- 
less, all data points lie on a line with a slope that is also con- 
sistent with the one we determined for the Msmbb — Mbulge 
relation in paper I. This indicates that the ratio of SMBH 
to halo mass does not depend significantly on the redshift. 

We also made an important assumption regarding the 
actual mass of the satellite systems in which the MBHs are 
embedded. The semi-analytical code we used to follow the 
satellite dynamics considers every satellite naked if its mass 
has dropped to something less than about 0.3 per cent of 
its original mass 3 . For all 'naked' MBH satellite systems 
arriving at the infall radius we have therefore assumed a 
mass equal to this limit. We do impose a lower limit of 3 x 
1O 3 M for the mass of any naked satellite system, roughly 
corresponding to 1 per cent of the mass of the minihaloes and 
large enough for the presence of a baryonic core of the order 
the MBH mass. It is this remnant material of the satellite 
associated with an MBH that allows MBHs to the delivered 



2 The SMBH to halo mass has also been determined on the basis 
of velocity dispersion data outside the bulge dominated regions 
IFerrarese 20021 

3 The actual criterion refers to a limiting tidal radius of a satel- 
lite system, beyond which it is considered 'naked'. This criterion 
roughly corresponds to the 0.3 per cent mass limit we used. 



to the centre efficiently by dynamical friction. This point 
has also been mentioned by Yu (2002) . 

For the halo of final mass 1.6 x 10 12 M we find the 
following. Of the satellite MBH systems crossing the infall 
radius about 96 per cent do arrive at the centre within a 
Hubble time. These MBHs contain 99 per cent of the mass 
of all MBHs in the infall region. If we reduce the satellite 
mass to 0.03 per cent of its original value and the lower limit 
to 3 x 10 3 M Q we still find that 88 per cent of the MBHs in the 
infall region (accounting for 97 per cent of the mass) end up 
at the centre within a Hubble time. In haloes of lower final 
mass virtually all MBHs systems crossing the infall radius 
travel to the centre within a Hubble time. 



2.3.2 Dynamical Friction and Gas Infall in Galaxy 
Mergers 

MBHs with masses and orbits that would not allow them 
to travel to the centre within a Hubble time may still end 
up in the centre if the environment changes in a way that 
allows dynamical friction to act more efficiently. Mergers 
between galaxies of similar mass ('major mergers') through 
the accompanying violent dynamical processes may induce 
the dramatic changes to the matter distribution in the 
central region that we are looking for. Hydrodynamical 
simulations of mergers of gas rich galaxies (Barnes 2002 
|Barnes fc Hernquist 1996|lNaab fc Burkert 200lTl show that 
up to 60 percent of the total gas mass of the two galaxies 
can end up within a region just a few hundred pc across at 
the centre of the merger remnant, thereby triggering star- 
bursts and initiating the fuelling of quasars 4 . This highly 
dense gaseous core should remain until feedback from star 
formation, supernovae or quasar activity begins to drive it 
out again. 

We have already identified gas infall as a potential 
mechanism for accelerating the evolution of a hard MBH 
binary. We now look at how this effect can also boost the 
efficiency of dynamical friction as MBHs move in the central 
region of the merger remnant. If the respective central- kpc 
region of each galaxy ends up within a kpc or so of the new 
centre of the remnant after merging, then so will any MBHs 
that were present within a kpc of the original galactic cen- 
tres. This seems reasonable since the central regions, sitting 
deepest in the gravitational potential of each galaxy, will be 
most stable to tidal disruption. It is plausible that MBHs 
originally present there travel along the line connecting the 
centres of the merging galaxies, towards the new global po- 
tential minimum of the merger remnant. 

This migration of MBHs could be accelerated by the 
rapid infall of large amounts of gas into a small core region. 
The question then is whether the density and the bulk in- 
flow velocity of the gas is sufficient to accelerate any orbiting 
MBH towards the centre via dynamical friction. The com- 
plex dynamical events and rapidly varying geometry make 
it difficult to determine accurately any radial bulk inflow 
velocities and gas densities as they evolve with time, and 
particularly so in the core central region, where things are 

4 The galaxy interaction leads to the shock heating of the gas 
and the formation of a bar in the gas distribution. Radiatively 
cooled gas then falls towards the centre along the bar. 
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made more difficult by the resolution limits of the simula- 
tions. Even a rough estimate of these quantities may already 
yield useful constraints for our problem, however. In the fol- 
lowing we try to estimate the timescale for infalling gas to 
carry two MBHs towards the centre by dynamical friction 
before the gas leaves the central region again. 

For two disk galaxies, each having a total mass M ga i = 
2.75 x 10 11 M Q and gas mass M gas = 5.5 x 10 9 M Q in the 
disk and initial pericentric separation of R p = 8kpc, Barnes 
& Hernquist (1996) find that roughly 60 percent of the total 
gas mass end up within a region approximately 100 pc in ra- 
dius some 1.5 Gyrs later. Assuming the centre of the merger 
remnant will be half-way between the galaxies at initial sep- 
aration and the gas flows towards this point along the line 
connecting the merging galaxy centres, the average inflow 
velocity of the gas is Vfi ow ~ 4kpc/1.5Gyr w 2.7 kms -1 . 
Further we assume that the gas density in the flow is of 
the same order as the final gas density in the centre, i.e. 
Pfiom ~7x 10 9 M o /(4 x 100 3 pc 3 ) w 1.6 x 10 3 Mpc pc -3 . 
Putting this into eg. 1171 we find that the MBHs get acceler- 
ated to within 1 per cent of Vfiow in a time 

T **~ ~ 2 x 106 {xtrk ) 1 yr (19) 

Even for seed mass MBHs, this time is negligible compared 
to the time scale of the inflow of about 1-2 Gyrs. Although 
we do not know what exactly will happen to the MBHs 
once they have entered the central 100 pc of the remnant, 
gas inflow seems at least capable of efficiently transporting 
MBHs from the outer radius of the infall region (~ lkpc, 
for galaxies of this size) down to at least 100 pc. This holds 
in particular for MBHs that previously would not have been 
able to spiral anywhere near the centre. 

2.3.3 Multiple MBH interactions 

Up to now we have implicitly assumed that a MBH binary 
merges before it forms a new binary with another incom- 
ing MBH. However, given the number of inspiralling MBHs 
involved it is possible that there are cases where incom- 
ing MBHs will find another MBH binary that has not yet 
merged. If the binary system and the incoming third MBH 
come close enough then there is a possibility that one of the 
MBHs will be ejected jBegelman, Blandford fc Rees 19801 - 
The required proximity of an incoming MBH and MBH bi- 
nary is facilitated primarily by major mergers. This will be 
the case particularly at high redshifts, as most MBH inter- 
actions will be between MBHs with similar masses (of the 
order of the MBH seed mass), which also implies similar 
masses of their associated haloes. 

As far as our model is concerned we have not incorpo- 
rated the possibility of triple interactions in our results. We 
assume that stellar and gas dynamical effects in the wake 
of major mergers always lead to swift merging of MBH bi- 
naries, primarily through accelerating evolution of the hard 
binary stage. Here we consider what happens if we drop this 
assumption. For the case that stellar and gas dynamics have 
no mitigating effect, we can determine a rough limit for the 
number of triple interactions and resulting sling-shot ejec- 
tions of MBHs. 

In figure|2|we show the number of mergers per unit time 
(10 6 yr) for model A and C. From top to bottom the curves 



are for 1.6 x lO 12 ' n ' lo M haloes. The dashed lines are the 
best fits. Let us consider the largest merger rate, which we 
find at z ~ 5 — 10 in the precursors of haloes with final 
mass 1.6 x 10 12 M Q . The rate is dJV mer9E / dt ~ 0.07 Myr" 1 
implying that two MBHs arrive at the centre within ~ 15 
Myrs. For the lowest rates AN merge / dt ~ 0.001 Myr -1 two 
MBHs would reach the centre within about a Gyr. These 
timescales would have to be compared with the hardening 
time of a MBH binary at the centre to determine the likeli- 
hood of a triple interaction. Yu (2002) finds that the hard- 
ening timescale can be significantly larger for equal mass 
binary systems, which again is more likely to be the case at 
very high redshifts. 

Although the central halo MBH merger rate as shown 
in figure |3] is highly uncertain at high redshifts, there is a 
clear declining trend towards high redshifts. Consequently 
we would expect multiple MBH interactions to be less im- 
portant. Overall, this is different from Volonteri et al. (2003) 
, who find that there is a very significant probability for 
triple interactions. One important reason for this discrep- 
ancy is that they explicitly did not take into account tidal 
stripping of satellite systems and the resulting increase in 
the dynamical friction infall timescale. When MBHs are in- 
volved in triple interactions, they find that those MBHs that 
get ejected MBH will actually leave the galaxy (remnant) al- 
together (Volo nteri, Haardt fc Madau 20 03). These ejected 
MBHs could thus constitute a population of inter-galactic 
MBHs. Due to the nature of the process these MBHs will be 
mostly seed mass MBHs and will have been stripped naked, 
so that they will be virtually impossible to detect in the 
IGM. 

The overall result is then, that while the formation of 
multiple MBH systems is possible, the probability of such 
events is likely less important than analytical arguments sug- 
gest. Previous estimates typically did not account for the 
tidal stripping of infalling satellites let alone their detailed 
dynamical evolution, but also the role of gas dynamics at the 
centre. The first effect implies that the rate at which MBHs 
arrive from the outer parts of the halo at the halo centre is 
lower, whereas the latter means that once MBHs have ar- 
rived at the centre, they form a binary that coalesces faster. 
In the light of this, we believe our assumption of efficient 
merging to be justified, particularly at redshifts <; 15 where 
MBHs mergers actually produce a strong enough gravita- 
tional wave signature that could be detected as we will see 
in the next section. 

2.4 Rate of SMBH-MBH mergers 

We now attempt to obtain the cumulative signal of all merg- 
ers. To do this we essentially need to integrate over the mass 
function of haloes times the number of mergers per halo and 
finally integrate again over the relevant redshift ranges. 

In what follows, we will only consider models A and C. 
In model D the rate of mergers is too low to derive a mean- 
ingful estimate of the average merger rate per halo, and 
consequently the total rate of events received from across 
the sky. Model B only differs from A in that it has a larger 
seed MBH mass. We do not expect this to make a signif- 
icant difference when it comes to mergers with the centre. 
When calculating the dynamical friction time from the point 
the MBH systems cross the infall radius (cf. section r2.3.H . 
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Figure 3. Rates of halo MBHs merging with central (S)MBH, as a function of redshift. Results are shown for models A and C and for 
final halo masses of 1.6 X 10 12 ' 11 ' 10 M Q (three sets of curves, from top to bottom). The dashed curves represent the same best fits to the 
data scaled by a factor corresponding to the different final halo mass (see text). 
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Figure 4. Rates of halo MBHs merging with central (S)MBH, as a function of redshift. Results are shown for models A and C and 
final halo mass of 1.6 X 10 12 Mq with the dashed line representing a best fit. From top to bottom, the three sets of curves show the 
contribution of events with dimensionless strain amplitudes above a given value as shown. 



we assumed the MBH satellite system's mass to be ~ 0.3 
per cent of the original satellite mass. The difference in the 
satellite's MBH mass in model A and C is thus insignificant. 
The larger MBH mass in model B will also only lead to very 
slightly increased strain amplitudes when MBHs eventually 
merge. 



The number of MBH mergers per halo is shown in figure 
03 for models A and C. The dashed lines represent a single 
cubic-polynomial fit in log space that is multiplied with a 



simple scaling factor for different halo masses. 



dre merge (Mfc,o, z) 
dt 



1.6 x 10 12 M Q/ 

b +b! z + b 2 z 2 +b 3 z 3 



Myr (20) 

where Mh,o denotes final halo mass. For model A we find 
(bo,bi,b 2 ,b 3 ,/3) = (-3.75,0.12,0.022,-2.63 x 10" 3 ,|), 
while for model C we find (-3.52,0.145, 7.5 x 10~ 3 ,-1.57x 
10 -3 , i). In this form the merger rate can now simply be 
divided into a product of a redshift dependent part, u(z), 
and another one scaling with final halo mass 



u(z) (21) 



dn„ 



dt 
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Figure^Jshows the merger rate per halo with strain am- 
plitudes h c > 10 -21 and 3 x 10 -21 for a final halo mass of 
1.6 x 10 12 Mq only. While the best fit parameters change ac- 
cordingly, we have assumed that the scaling with halo mass, 
as parametrised by /3, remains the same. Although we have 
not explicitly stated any error bars, uncertainties will be 
quite large, since the merger rate per halo is determined by 
averaging over individual merger events occurring over large 
timescales. The only way of reducing this uncertainty is to 
run a significantly larger number of merger trees. 

When determining h c , we assume that the merging 
MBHs emit waves at the maximum frequency which is only 
expected to occur very briefly towards the end of the inspi- 
ral phase due to the rapid sweeping up in frequency. This, 
however, should not affect our result significantly, since h c 
depends only weakly on the frequency (cf. equation 10). 

Using the merger rate per halo we can already deter- 
mine as a function of epoch a lower limit to the rate, dv r ec 
of received SMBH-MBH merger events. 

dv rec (z) I 



dV 



> 



1 


+ z 




poo 


X 






Ju h 


u 


(*) 


1 


+ z 






X 






J M h 



dnh dn„ 



dMfc.o 



dn h 



dt 



-(M h , ,z) dM h , 



Mi, 



dM h ,o V 1.6 x 10 12 M G 



dM, 



(22) 

where V denotes volume and (1 + z) accounts for the 
redshift induced reduction in the rate of received events. 
dnh/dMhfi is the differential mass function of haloes of mass 
M hi0 



dn h 
dM hx 



2 Po 

7T M h c 



diner 



<?(M h fi 



-exp 



2a(M h 



d lnM hi0 

where po is the present-day cosmic matter density, a(Mh,o) 
is the variance in the linear matter density field on scales 
corresponding to mass Mh.fi, and 5 C is the critical overden- 
sity for collapse. 

Equation l22l is a lower limit, as it only considers mergers 
in the main trunk of the tree that grows to become a single 
halo of mass Mh at z = 0. To account for the mergers in side 
branches before they become incorporated in the main trunk 
of the tree the Press-Schechter (PS) mass function, which 
we had only computed for final halo masses above, needs 
to be replaced by the redshift dependent mass function. In 
addition we need to re-express the redshift dependent MBH 
merger rate per halo in terms of the halo mass at the redshift 
concerned. From our simulations we found that the average 
mass of a halo at redshift z, given that it grows into a halo 
of final mass Mh,o today, can be approximated by 

M h (z) = M h (M h<0 ,z) « exp [12 + m - 0.04 rh z] (23) 



where rh = 11.5 + In 



We can now substitute 



Mho 
s 1.6xl0 lu M Q / 

this into equation 1211 to eliminate the dependence on Mh,o 
and obtain 



dn„ 



dt 



-(M h ,z) = 



M h 



1.6 x 10 5 M 



P 



u(z) 
(24) 



The lower integration limit in ea 1221 is determined by 
the minimum mass of haloes that harbour central MBHs. 
This mass has to be at least of order 10 7 Mq, since we 
argued above that only at this mass has every halo acquired 
one seed MBH. This is important since we have implicitly 
identified the number density of haloes with that of central 
MBHs. 



If we now multiply 



div 



with the volume of a spher- 



ical shell at radius D(z) we obtain the total rate of events, d 
Urec(z), received from sources at distance D(z) to D(z)+dD 

dv rec (z) = AVr *y {z) ^D{zfdD (25) 

The resulting differential redshift distribution of the received 
rate of events is shown in figurc|S] where we have also shown 
the contribution from events with associated dimensionless 
strain amplitudes h c above 10 -21 and 3 x 10 -21 . The latter 
is based on the fits to the merger rate per halo shown in 
figure^] Again, the errors are quite large and are difficult to 
quantify precisely, primarily because of the uncertainties in 
the merger rate per halo to which we now also have to add 
the uncertainty in the halo mass scaling, which we described 
by a power law with index /3. 

This differential distribution can then be integrated over 
all redshifts to obtain the total number of merger events 
received. 



2.5 Detections and the distribution of strain 
amplitudes 

For actual detections with LISA, we need to take into ac- 
count its restrictions on the minimum detectable strain am- 
plitudes. For the MBH masses under consideration here, and 
the resulting wave frequencies, LISA should be able to see 
events with h c ~ 10 -21 , and may even see events as weak 

0.01 Hz. For the combination of 



as h c ~ 10~ near f g 
MBH masses and the redshift range considered, all events 
yield strain amplitudes larger than 10 -23 , and LISA may 
therefore be able to detect most if not all of these events. 
Adopting the more conservative LISA limit of h c ~ 10~ 21 
the number of detections is somewhat smaller as shown by 
the corresponding curve in figure [5] 

We are interested in detections of events at cosmolog- 
ical distances, and to do so the conservative LISA detec- 
tion limit requires a minimum chirp mass of about 10 4 Mq 
IIHaehnelt 199411 . MBHs of this size, however, will typically 
be hosted inside haloes of at least 10 9 Mq, which could 
therefore be taken as the effective lower limit of integration 
in eauation l22l In fact dv rac is not very sensitive to the par- 
ticular choice of integration limits. The event rate of mergers 
in haloes as determined by our simulations is very small at 
the low-mass end, which therefore does not contribute sig- 
nificantly to the overall rate. Similarly, the steep decline in 
the mass function at high masses imposes an effective upper 
limit. 

The redshift distribution of merger events can now be 
integrated across the redshift range to obtain the total num- 
ber of events received. The result is shown in figure^ where 
we have shown the total number of events above a given 
strain amplitude. As in the previous section, we have used 
a fit to the merger rate per halo above a given strain ampli- 
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Figure 5. Redshift dependent rate of received merger events for models A and C. From top to bottom, the three sets show the contribution 
from all events and those with dimensionless strain amplitudes above 10 — 21 and 3 X 10 — 21 respectively. 



tude to determine the merger rate per unit redshift. From 
this we see that LISA would be expected to detect some 
10 4 - 10 5 events per year. The difference between model A 
and C is too small to be detectable in practice, as the overall 
uncertainties in our results are of at least the same order. 

We have seen above that the strain amplitude is rather 
insensitive to the frequency of the waves emitted. The fre- 
quency does matter for LISA detections as its peak sensitiv- 
ity lies in a rather narrow band around 0.01 Hz. Merging seed 
MBH binary systems have very low strain amplitudes, es- 
pecially at high redshifts, but peak frequencies well above 1 
Hz. The only way of detecting these is during a significantly 
earlier phase of their coalescence when their frequency is 
lower. 



3 GAMMA RAYS FROM MBH FORMATION 
AND MERGING 

Another exciting possibility is that the more massive of the 
first stars forming at high redshift conclude their lives in 
hypernovae, which be visible as gamma-ray bursts (GRBs). 
GRBs should be detectable to high redshift by upcoming 
satellite missions, which would allow us to study events from 
the first era of star formation directly. In contrast to the 



other observations mentioned above, such observations can- 
not be used to verify the hierarchical MBH merging scenario. 
GRB observations might provide constraints on the abun- 
dance of the first massive stars, however, and in this way 
allow us to test this part of our model directly. 

GRBs are the brightest explosive events in the Universe. 
Typically they involve the emission of large amounts of ra- 
diative energy, of the order ~ iq 52-54 er g S m j rays, over a 
period of only a few seconds. A GRB is followed by an af- 
terglow in the X-ray, optical or radio parts of the spectrum 
and lasts from days to weeks. The total energy emitted in 
an afterglow is typically two or three orders of magnitude 
below that of the actual burst. Spectroscopy has been car- 
ried out on the afterglow emissions as well as the emission 
of what appear to be host galaxies of the GRBs, once the 
afterglow itself has faded away sufficiently. In this way it has 
been possible to constrain the redshift of a number of GRBs 
below z < 3. If GRBs could be used as standard candles in a 
similar way as type la supernovae, this would facilitate their 
use to calibrate the GRB distance-luminosity relation. On 
this basis one could then determine the maximum redshifts 
at which GRBs could still be detected by present and future 
observing missions, such as with the BATSE detector on the 
Compton Gamma-Ray Observatory or the SWIFT satellite 
due to be launched late in 2003 ( Lamb & Reichart 2000 ) . 
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Figure 6. Distribution of dimcnsionless strain amplitudes for events received from across all redshifts. 



For a number of GRBs for which redshifts could be estab- 
lished, Lamb & Reichart (2000) found that SWIFT could 
have in principle detected these out to redshifts of z ~ 20, 
and some to z ~ 70. For further details of GRBs their after- 
glows the reader is referred to the review by Piran (2000) and 
the extensive work by Bloom (2002) and references therein. 
For the following we assume that GRBs are in principle de- 
tectable out to redshifts beyond z > 20. 

A range of processes have been discussed that could be 
the cause for GRBs. As mentioned, one of these is the pos- 
sibility that GRBs are produced by the collapse of massive 
stars. In the following we exclusively focus on this possibility 
and attempt to estimate the number of GRBs from massive 
population III stars in our model. 



3.1 Number of GRBs from massive star collapse 

The procedure to compute the total number of GRBs across 
all relevant redshifts is very similar to the one we used to de- 
termine the gravitational wave events. It is different in that 
we only deal with the very first stars, and are not concerned 
with the subsequent merging of the MBHs these stars leave 
behind. This is therefore a purely analytical argument with 



no recourse to our simulation data concerning the merging 
and dynamical evolution of MBHs. 

Again, we estimate the abundance of minihaloes at the 
redshifts in question (z ~ 25). Depending on how many 
massive stars form per halo, we can then determine the ex- 
pected number of GRBs. A similar route has been followed 
in previous work. In this case the abundance of massive 
stars per halo/galaxy was assumed to trace the star forma- 
tion rate at the respective redshift IBromm fc Loeb 20021 
|Roy Choudhury fc Srianand 2002| . We assume that every 
minihalo between a maximum redshift z maXl corresponding 
to a 2.5-3.5cr peak, and a minimum redshift determined by 
the end of first star formation, which we take as z m in ~ 20, 
produces one star that explodes as a hypernova with an as- 
sociated GRB. 

We have seen that a hypernova leads to the expulsion 
of all gas from minihaloes. Each minihalo can thus only 
ever produce one GRB. This is a problem as we now re- 
quire each halo to retain a memory of whether and when 
it produced a GRB. For this reason it is not strictly cor- 
rect to use the Press-Schechter mass function to compute 
the number density of haloes at a given redshift and equate 
the result with the number density of GRBs. A halo more 
massive than M m i n may have already lost all its gas through 
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a hypernova at some earlier stage when it was significantly 
less massive. In practice, however, the Press-Schechter ap- 
proach is more than adequate for two reasons. First, the 
PS mass function is steepest for large masses correspond- 
ing to haloes collapsing from high peaks as is the case for 
all minihaloes in our redshift range. The number density 
of haloes with mass Mh ^ M m in is therefore vastly domi- 
nated by haloes with masses very close to M m in- Secondly, 
the survival time 5 for haloes collapsing from high peaks 
is very short QLacey fc Cole 1993| l. That means that by far 
most haloes with Mh ^ M m i„ will have formed at or very 
shortly before the redshift in question and the probability 
of previous star formation and a GRB is very small. 

If we identify a GRB event with the collapse of a mini- 
halo, the rate of GRBs at any given redshift is then approx- 
imately 



' i;V;/,B dz^A7vD(z) 2 dD(z) 



dz 



d 2 n h 
dt dM 



dM (26) 



where we have used the rate of change of the minihalo mass 
function 

d 2 n h 



dt dM h 
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d InMh 
1 d« 
u dt 
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~dl 



(27) 



dnh du 
dM h ~~dt 



and u = 5 c /[a(Mh) * D grom (t)]. We have shown the mass 
function and its rate of change in figure 

3.2 Rate of GRBs received 

For GRBs at the very high redshifts we are considering, 
no afterglow emission can be detected anymore. An ac- 
curate redshift determination is thus not possible unless 
a unique distance-luminosity relation for GRBs could re- 
ally be established. If all GRBs do have the same to- 
tal luminosity and peak photon flux, we have a unique 
distance-luminosity relation and thus a unique redshift- 
luminosity relation. Eqn. I27H above can then be con- 
verted into a flux distribution Fgrb(z) using dz oc 
-2tt j D( Z ) 3 /(10 52 - 54 erg)( dD/ dz)" 1 dF GRB . 

However, the case for a unique GRB distance- 
luminosity relation has yet to be made. To test our model 
we therefore focus on a prediction of the rate of GRBs to be 
compared with those detected by e.g. the SWIFT mission. 
We integrate eq. 1)27^ over the redshift range of population 
III star formation to arrive at the total rate of GRBs re- 
ceived. 



vgrb,. 



' {1 + z) -idi^B dz m 



This is shown in figure |H| together with rate of GRB 
events originating from any given redshift interval for three 
different values of M m in, assumed constant over the redshift 
range 20 < z < 30. We find that the total number of GRB 

5 The survival time is the time by which the mass of a halo has 
more than doubled. 




Figure 8. Number of GRB events received per year originating 
from redshift z to z + dz. The dashed curve shows the total 
number of GRB events from z 20. From top to bottom the 
curves correspond to a M m i n = {1,3, 10} X 10 5 Mq/i —1 . 



,„,„ = 10° to 10° 
5 M fl at z = 25, 



events received per year from this redshift range varies be- 
tween vgrb « 3 x 10 6 - 2 x lOV -1 fo- ™ - ln5 
Me respectively. Even for M m in — 3 x 10" IVl© at z 
it may still make more sense to consider a larger halo mass, 
since by definition baryons - and therefore stars - in mini- 
haloes only cool on a Hubble timescale, whereas we have 
assumed that GRB events coincide with the collapse of the 
haloes within which they occur. At z — 25, for instance, a 
halo mass three times larger than the minihalo mass we used 
would result in a cooling time that is about a tenth of the 
Hubble time at that redshift and so more in line with this 
assumption. 

A conservative estimate for the total number of GRB 
events received is vgrb ~ 2 x 10 5 yr _1 , assuming that halo 
collapse, star formation and GRB coincide and that first star 
formation ends at z = 20. How many of these could actually 
be observed depends on whether or not GRBs are luminous 
enough to still be detected at these distances. 



3.3 Gamma rays from MBH mergers 

If a small fraction of the gravitational energy from MBH 
mergers is released in the form of gamma rays - maybe due 
to the fast accretion of large amounts of material in the very 
close vicinity of the MBHs - we may be able to see these. 
For an inspiralling MBH binary the energy radiated away in 
gravitational waves, L aw is flShapiro fc Teukolsky 19851 
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Figure 7. Top: Redshift dependent comoving number density of haloes with mass > M m i n . From top to bottom the curves 
correspond to a M m i n = {1, 3, 10} X 10 Mq/i~ . Bottom: Rate of change of minihalo number density. The curves are in the same order. 



and combining this with the characteristic merger time scale 
eq|5]this becomes 

If we now multiply this with 1/fmax (eq|7J we can estimate 
the energy radiated away on the last orbit before final coa- 
lescence. 

^-"^ Garter- 

It is interesting to compare this with GRBs. A typical GRB 
lasts between a fraction of a second to minutes and emits to- 
tal energies between 10 51 - 10 54 erg l|Piran 2 000 ) and about 
a factor 100 less if we assume beamed emission. As far as the 
duration of GRBs is concerned this is matched by the period 
of the last orbits of MBH binaries with mass M12 ~ 1000 
M (r o = 1/fma* ~ 0.1 seconds) to M 12 ~ 10 6 M (r o ~ 100 
seconds). This means that if some fraction e 7 J> 3 x 10~ 2 of 
gravitational wave energy is emitted in gamma rays, merg- 
ing MBHs could account for some beamed GRBs. At the 
lower end the detection limit of the BATSE instrument on 
the Compton Gamma Ray Observatory 6 implies e 7 > 10 -6 
for MBH mergers to be detected out to cosmological dis- 
tances. For isotropically emitting GRBs e 7 would need to 
be a factor 100 or so higher. 



4 SUMMARY AND CONCLUSIONS 

In this paper we were mainly concerned with probing our 
model at high redshifts. The merging of MBHs in the con- 
text of the hierarchical merging of galaxies and haloes would 
lead to the emission of gravitational waves, that could be de- 
tected with the LISA mission due to be launched in the next 
decade. We can draw the following conclusions: 

B BATSE's burst sensitivity is quoted as 3 X 10 — 8 erg cm" 2 . 



• The merger rate per unit redshift peaks at z ~ 4-6. 
The peak shifts to lower redshifts if events with higher di- 
mensionless strain amplitudes are considered. 

• The predicted merger events all have strain amplitudes 
above h c > 10 -23 , and most above 10 -22 , so that LISA 
should be able to detect most of these. 

• The total number of MBH merger events received are 
about 10 4 -10 5 per year. This is at least an order of magni- 
tude larger than typical estimates for ('major') mergers of 
SMBHs in hierarchical structure formation scenarios. Our 
larger number arises from the inclusion of a lot of mergers 
with fairly large binary mass ratios, i.e. what we may qual- 
ify as 'minor' mergers. The required population of MBHs 
are the (merged) remnants of the first massive stars in the 
universe. An estimate of this population was the subject of 
paper I in this series. 

• In contrast to previous analytical work our predictions 
are based on a more accurate and explicit modelling of the 
dynamical evolution of MBHs in galactic haloes. This has al- 
lowed us to follow MBH merging in haloes at higher redshift 
and therefore also down to lower precursor masses. 

• If the collapse of the first massive population III stars 
into MBHs is accompanied by gamma-ray bursts (GRBs), 
these might allow us to probe the abundance of the first 
stars as well as their spatial distribution. We predict 10 - 
10 6 GRBs per year originating at redshifts higher than 20, 
if there is one massive star per star-forming minihalo asso- 
ciated with a GRB. 

• A fraction of the gravitational wave energy in MBH 
mergers may be released in gamma rays, e.g. through in- 
teraction with remnant gas in the close vicinity of the co- 
alescing binary. If GRB emissions are beamed this fraction 
must be e 7 > 10 -6 for the resulting bursts of gamma rays to 
be detectable at cosmological distances and e 7 > 3 x 10 -2 
to possibly account for some of the proper GRBs. However, 
the short event duration makes it difficult to identify them 
unless a GRB style after-glow is associated with them. 

While the prediction of gravitational wave events and 
GRBs both depend sensitively on the initial assumption of 
massive star formation, the additional assumptions of effi- 
cient merging for gravitational waves and that GRBs are ere- 



MBH remnants of the first stars III: observational signatures from redshift evolution 13 



ated in the process of massive star collapse, further increase 
uncertainty. All these assumptions are in principle verifiable 
primarily with improved numerical simulations (The role of 
gas in BH binary mergers, improved resolution in star for- 
mation and collapse simulations in a cosmological context). 
For the time being, the prediction of event rates of 10 to 
100 a day should provide motivation to look out for them in 
present and upcoming observations. Qualitatively the most 
important thing is to note, that we could potentially observe 
the highly energetic events accompanying both the forma- 
tion (via GRBs) and merger (by gravity waves) of MBHs, 
and therefore probe the high redshift regime. 
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